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We examine numerically and qualitatively the LemaTtre-Tolman-Bondi (LTB) inhomogeneous 
dust solutions as a 3-dimensional dynamical system characterized by six critical points. One of the 
coordinates of the phase space is an average density parameter, (£2) , which behaves as the ordinary 
Q in Friedman-Lemaitre-Robertson-Walker (FLRW) dust spacetimes. The other two coordinates, a 
shear parameter and a density contrast function, convey the effects of inhomogeneity. As long as shell 
crossing singularities are absent, this phase space is bounded or it can be trivially compactified. This 
space contains several invariant subspaces which define relevant particular cases, such as: "parabolic" 
evolution, FLRW dust and the Schwarzschild-Kruskal vacuum limit. We examine in detail the phase 
space evolution of several dust configurations: a low density void formation scenario, high density 
re-collapsing universes with open, closed and wormhole topologies, a structure formation scenario 
with a black hole surrounded by an expanding background, and the Schwarzschild-Kruskal vacuum 
case. Solution curves (except regular centers) start expanding from a past attractor (source) in the 
plane (Q) = 1, associated with self similar regime at an initial singularity. Depending on the initial 
conditions and specific configurations, the curves approach several saddle points as they evolve 
between this past attractor and other two possible future attractors: perpetually expanding curves 
terminate at a line of sinks at (£2) = 0, while collapsing curves reach maximal expansion as (fi) 
diverges and end up in sink that coincides with the past attractor and is also associated with self 
similar behavior. 

PACS numbers: 12.60.Jv, 14.80.Ly, 95.30.Cq, 95.30.Tg, 95.35.+d, 98.35. Gi 



I. INTRODUCTION 

Inhomogeneous dust solutions with spherical symme- 
try are among the oldest, simplest and most useful ex- 
act solutions of Einstein's equations. These solutions 
were initially derived independently by Lemaitre (1933) 
and Tolman (1934) and then re-derived by Bondi (1947), 
hence they are known as the Lemaitre-Tolman -Bondi 
(LTB) solutions (see [I] and references quoted therein 
for a comprehensive review on these solutions and their 
applications). 

In practically all applications of these solutions the 
standard original variables in which they were derived 
are used, though alternative variables amenable to an 
"initial value" treatment have been proposed (see [5] and 
references quoted therein). While not the only type of 
parametrization for these solutions, we find these vari- 
ables particularly useful for a numerical treatment. 

We propose in this paper to examine these models 
within the framework of numeric and qualitative tech- 
niques known generically as "dynamical systems" 
This approach requires as a first necessary step express- 
ing the evolution equations as a proper system of first or- 
der autonomous ODE's (ordinary differential equations). 
This necessary requirement would, apparently, exclude 
inhomogeneous LTB solutions because their evolution 
equations are necessarily PDE's (partial differential equa- 
tions). However, under a "3+1" covariant decomposition 
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based on the 4- velocity field [I], the evolution equations 
for LTB dust solutions become first order autonomous 
PDE's containing only time derivatives [3]. We prove 
in this case how the spacelike gradient equations (con- 
straints) are automatically satisfied for all times once 
initial conditions are selected so that these constraints 
hold in an arbitrary initial hypersurface 3 7^ of constant 
time. Rigorously (see Appendix B) , the solution curves of 
these evolution equations are then equivalent to a proper 
dynamical system built with ODE's, under the strong re- 
striction of complying with the special set of initial con- 
ditions that is compatible with the spacelike constraints 
of the PDE system. 

In order to describe LTB dust solutions as a dynamical 
system we re-write the first order evolution equations 
of the 3+1 decomposition [1] in terms of initial value 
variables defined in reference [2]. These variables are 
introduced and generalized in section III and in sections 
IV and V we show how they lead in a very natural way to 
a system of three autonomous evolution equations that 
is similar to those obtained with "expansion normalized" 
variables for other known space-times (FLRW, Bianchi 
models, Kantowsky-Sachs, see [3])- 

Initial conditions (section VI) are selected from two 
"primitive" initial value functions (density, scalar 3- 
curvature) defined along an initial and regular hyper- 
surface, 3 7l marked by constant comoving time. These 
functions determine the initial values for the phase space 
variables in such a way that the two spacelikc constraints 
characterizing the 3+1 evolution equations for the LTB 
solutions are satisfied at the z % and at subsequent 3 T. A 
third initial function determines the number of Regular 
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Symmetry Centers (RSC) at the 3 % (see Appendix A), 
and thus it determines the topological class of the 3 T, 
leading to hypersurfaces with "open" topology, (home- 
omorphic to M 3 , one RSC); "closed" topology (homeo- 
morphic to S 3 , two RSC) "wormhole" (homeomorphic to 
§ 2 x K or to § 2 x S 1 , zero RSC). 

The resulting 3-dimensional phase space is discussed 
in section VII, with its critical points, invariant subsets 
and particular solutions given in section VIII. This space 
is parametrized by three functions: a dimensionless shear 
scalar, S, a density average parameter, (p.), and a den- 
sity contrast function A^ m K As long as initial condi- 
tions are selected so that unphysical shell crossing sin- 
gularities do not arise [21 [7] , these variables are either 
bounded or can be trivially compactified for all configu- 
rations. These functions have a straightforward physical 
interpretation: the average (ft) behaves exactly as the 
standard ft in a FLRW dust spacetime: its initial val- 
ues (ftj) determine the dynamical evolution of each so- 
lution curve: re-collapse (if (ftj) — 1 > 0) or perpetual 
expansion (if (ftj) — 1 < 0). Thus, (ft) defines invari- 
ant subspaces (ft) = 0, (ft) = 1, < (ft) < 1 and 
(O) > 1. The density contrast function A^ m ) provides a 
non-local characterization of the type of inhomogcneity 
along the 3 T (rest frames of fundamental observers): a 
"clump" if -1 < A( m ) < or a "void" if A^" 1 ) > 0, 
with A^" 1 ) = at the RSC. This variable also defines an 
invariant subset: A' 1 ™) = — 1, corresponding to vacuum 
Schwarzschild-Kruskal solutions. Another invariant sub- 
space is given by the line S = — marking the 
homogeneous and isotropic FLRW subcase, containing 
as subcases the Einstcin-de Sitter and Minkowski-Milne 
universes. 

We examine in detail in section IX the phase space 
evolution of several representative dust configurations: 
low density perpetually expanding ( "hyperbolic" dynam- 
ics) with (ftj) — 1 < 0, as well as a case of void forma- 
tion scenario in which an initial density clump becomes 
a void (a numeric realization of configurations proposed 
in [10]). We also examine configurations with high den- 
sity re-collapsing ("elliptic") dynamics with open, closed 
and wormhole topologies ((ftj) — 1 > 0); a "structure 
formation" scenario: solution curves near the RSC re- 
collapse into a black hole (ftj) — 1 > while "exter- 
nal" curves perpetually expand into a cosmic background 
(ftj) — 1 < 0. We also examine the vacuum subcase, which 
is the Schwarzschild-Kruskal spacetime given in comov- 
ing non-static coordinates made up by radial timelike 
geodesies. 

All solution curves (except the RSC) for all configura- 
tions start their evolution at a past attractor (source) at 
(ft) = 1, S = 1/2, A(" 1 ) = -1 that represents an initial 
big bang like singularity and is associated with a self sim- 
ilar regime. The evolution of each solution curve depends 
on the sign of (ftj) — 1 and is fully contained in the invari- 
ant subsets associated with (ft). Curves with (ftj) — 1 = 
( "parabolic" evolution) remain in the plane (ft) = 1 and 
terminate in a future attractor (sink) associated with the 



Einstein-de Sitter zero spacial curvature FLRW universe. 
Curves with (ftj) — 1 < ("hyperbolic" evolution) evolve 
towards a future attractor given by a line of sinks at 
(ft) = S = 0, while curves with (ft,) - 1 > ("elliptic" 
re-collapsing evolution) evolve into a maximal expansion 
state (ft) — * co and then terminate at the same critical 
point from which they started, though this point is now 
a sink or a future attractor. 

The RSC's of the configurations evolve separately from 
the rest of the curves along a line with A( m ) = S = 0, 
starting always from a source at (ft) = 1 and going to- 
wards a sink at (ft) = or towards (ft) — * oo, respec- 
tively, for expanding and re-collapsing configurations. In 
all cases some of the curves pass near a saddle close to 
the RSC associated with homogeneity (S = A^ m ' = 0), 
while in the structure formation scenario some curves 
also approach another saddle point that splits collapsing 
curves (ftj) > 1 from expanding ones (ft;) < 1. The solu- 
tion curves of the vacuum case are confined to the plane 
A( m ) — — 1, evolving towards (ft) = or (ft) — > oo, de- 
pending on whether the radial geodesies are bound or 
not. 

We summarize and discuss the results from the article 
in section X. We also provide three appendices dealing 
with important issues: RSC's and geometric properties 
of hypersurfaces 3 T (Appendix A), treatment of evolu- 
tion equations given as PDE as a restricted dynamical 
system (Appendix B) and analytic solutions in terms of 
the variables used in the article (Appendix C). 



II. LTB DUST MODELS IN THEIR ORIGINAL 
VARIABLES. 

The Lemaitre-Tolman-Bondi (LTB) metric [HH1E1E] 
is the spherically symmetric line element 

ds 2 = -c 2 dt 2 + dr 2 + Y 2 (d8 2 + sin 2 9dcj) 2 ), (1) 

1 — K 

where Y = Y(t,r), Y' = dY/dr and K = K(r). The 
momentum-energy tensor usually associated with (JlJ is 
that of a dust source: 

T ab = pc 2 u a u\ (2) 

where p = p(t,r) is the rest-mass density (in gm cm -3 ). 

For the metric ([I]) and source ^ with a comoving 4- 
velocity u a — 6ft, Einstein's field equations G a b = kT°" 
with k = 87rG/c 4 , yield: 



2M' = npc 2 Y 2 Y', (4) 

where M = M(r) has units of length and Y = dY/dx°. 

The only nonzero kinematic parameters are the ex- 
pansion scalar G = u a - a and the shear tensor: a a i, — 
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V( a; b)U - (Q/3)h a b, given by 

2Y Y' 



6 = 



Y 



a a b = diag[0,-2£,E,£], 



Y 
Y 



Y' 



(5) 



(6) 



Other important quantities are the "Electric" Weyl ten- 
sor: E a i, — C a bcd,u c u d and the 3-dimensional Ricci ten- 
sor, 3 1Z, of the hypersurfaces orthogonal to u a (marked 
by constant values of x° = ct): 



E a b = diag[0, 

M' 



-2S,£,£], 
M K 



3Y 2 Y' Y 3 
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pc 



M 



2(KY)' = 3 KY 2 Y'. 



(7) 



(8) 



The usual treatment of LTB dust solutions consists in 
solving analytically the evolution equation ([3]). The so- 
lutions are usually classified in terms of the sign of K, 
which determines the type of evolution of dust layers: 
perpetually expanding/collapsing "parabolic" (K = 0) 
and "hyperbolic" (K < 0), and re-collapsing "elliptic" 
(K > 0). These analytic solutions are given in Appendix 
C. However, instead of using these solutions, we will fol- 
low in this article a different approach based on defining 
new variables more suitable for a qualitative and numer- 
ical analysis. 



III. VOLUME AVERAGES AND SCALING 
LAWS 

The form of p and 3 TZ in terms of 2M' and (KY) ' 
in Q and ^ suggests considering the following volume 
averages along an arbitrary hypersurface 3 7^ orthogonal 
to u a and marked by constant x\ = cti 



{A) 



JdVi 



dV l = Y 2 Y( sin 0drd9d(j), (9) 



where the subindex j denotes evaluation at t = tj and we 
have assumed that 3 7l is fully regular in the integration 
range. Applying ^ to Q and |8]) we get 



KC 

~3~ 
1 



(P) 



KC 



l n) = 



J pY 2 Y'dr _ 
J Y 2 Y'dr ~ 
1 J 3 1lY 2 Y'dr 
6 J Y 2 Y'dr 



2M 
K 

- y2- 



(10) 



(ii) 



where we have dropped the subindex j and have taken 
the lower bound of the integration range of r to be deter- 
mined by a suitable boundary condition on pi and 3 lZi, 
for example a RSC. See Appendix A and section IX-E. 



The averages (p) and ( 3 1Z) are non-local quantities de- 
pending on the form of p and 3 1Z along the integration 
range. We define the following "contrast functions" com- 
paring these quantities with their local counterparts: 



P-(P) 
(P) ' 

3-7? _ IM 



w = 3 n-( 3 n) 

( 3 TZ) ' 



p=(p)[l + A("% (12) 

> 3 K= { 3 K)[1 + A^}, 

(13) 

The interpretation of these these contrast functions [5] 
follows by using the definitions Q, (10 1 and (12 1 and 
integrating by parts along the 3 T, leading to 



KC 

6M 



J p'Y 3 dr, 



(14) 



where this integral is evaluated from a RSC. Since M > 
and Y > and assuming p > 0, the only quantity that 
can change sign is p' , therefore, looking at the radial 
variation of a density profile from the symmetry center 
along an arbitrary 3 T we have 



<0 o p'<0, p<(p), 
>0 ^ p'>0, p>{p), 



density clump 
density void 



(15) 



Likewise, for the contrast function A^ k \ using the def- 
initions Q, (111 and ( 13 ) and integrating by parts along 
the 3 T we obtain an analogous expression to (14): 



A« =— / 3 n'Y 3 dr 1 
6it J 



(16) 



However, while a negative p is not physically interesting, 
there is no physical objection for 3 1Z being positive or 
negative, or changing sign in the allowed range of r. If 
3 1Z > for all the allowed range of r, we have the same 
results as with A (m) in (151, but if 3 TZ < 0, then curva- 



ture clumps or voids are defined by the opposite signs [5] . 

We define now the following "scale factors" related to 
the metric functions: 



Y 

Y 1 

Y'/Y 

Y[jYi 



£'/£ 
Y(/Yi 



(17) 
(18) 



relating Y and Y' evaluated at an arbitr ary and a fidu cial 
(or "initial") hypersurface 3 %. Since (10 1 and (111 are 



valid at any T, we obtain the scaling laws 



(P) 
{ 3 K) 



£3 ' 

£ 2 



while from Q, ^ and ([l7|-([l8]| we get 

Pi 

P £3 F ' 

1 1 



! ^)(i-r) 



(19) 
(20) 

(21) 
(22) 
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Comparing (12 



get the scaling 



I-([13! with @-(|20} and (|2T]h([22| we 

aws for A( m ) and AW 



1 + A (m) = 



1 + A 



(to) 



r 



2 (fe) = 2/3 + A W 

3 r 



(23) 
(24) 



The quantities 0, £ and £ given in ([5]), Q and 
take the forms 



e = 

£ = 



u r 
T + f : 
r 

~3T' 



(25) 
(26) 



f = y[p-(p)] = T WAW, (27) 



the evolution equation ([3]) becomes 



(28) 



while the LTB metric (|T|) takes the Friedmanian form 



ds 2 = 



l di 2 



T 2 Yp dr 2 

r 



17 



sin 2 ( 



(29) 



(p)' 


3F' 


a' 1 ™- 1 = 




(p> 








3y 


A (fe) = 











It is also helpful to express the spacial (i.e. radial) gra- 
dients of (p) and ( 3 1Z) in terms of the contrast functions 
and T. With the help of @ , ([§ , ([lOj) , ([II]), ([12} and ([13} 
we obtain the following useful relations: 



(30) 
(31) 



Notice that these equations are strictly valid along all 
hypersurfaces 3 T. 



IV. 3+1 DECOMPOSITION. 

From a covariant 3+1 decomposition of spacctimc 
based on the 4-velocity field u a presented in [31 H] , the 
field and energy balance equations for a dust model char- 
acterized by (1} and ([2} are equivalent to the following 
set of first order evolution equations: 



together with the constraints: 



2 ~ 



V fc a" a - -v„ e = 0, 

KC 2 

V b E b a - — V a p=0, 



(33a) 
(33b) 



where A = u c V C A is the convective derivative along 
the 4-velocity, V a A bc = h a d VdA bc is the spacelike gra- 
dient, tangent to the 3 T hypersurfaces and orthogonal 
to u a , and A^ = A^ - (l/3)A c c h ab is the spacelike, 
symmetric, trace-free part of a tensor A ab . 



Equations ( 32a I— ( 32d I are 4 tensorial equations for 4 
tensorial quantities O, p, a a bl E a b , which together with 
the constraints ( 33a I— ( 33b I form a completely deter- 



mined system of partial differential equations (equiva- 
lent to Einstein's field equations). Bearing in mind that 
u a — 8 a o and considering the forms of the trace-free ten- 
sors <7 a b and E a b in ([6} and Q, the system (32 1 reduces 
to the following set of scalar equations: 



e 

3 




2 

KC 

-^p = o, 


(32a) 


e 

3" 


& (ab) = 


o 

— Oa ab - (j l - a c a h)c 
3 


- E ab = 0, 


(32b) 


£ = 


fi{ab) _ 


-^pa ab -QE ab H 


-3CT (a c £ b)c , 


(32c) 


£ = 


P = 


-p&, 




(32d) 


P = 



e 

3 

2 

~3 



2E 2 - 

es + s 2 -£, 

K 



KC 



pc l + 3£ 



-pQ, 



(34a) 

(34b) 

(34c) 
(34d) 
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and ( 33a H 33b I become 

T 



E ' 



Y' 

3E- = 0, 



/ 2 
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3£ 



Y' 
Y 



0. 



(35a) 
(35b) 



It is straightforward to express the system ( 34 I and the 



constraints (351 in therms of the new variables introduced 



in the previous sections. Considering (121, ( 23 ) and ( 25 I 
( 27 I and and defining 



i 3 



(36) 



the system (|34|) becomes the following set of equivalent 

(37a) 



evolution equations: 

H = -H 2 



(i 



(P), 



E = 3E 2 -2#E + 
(/>>• = -3 (p) H, 
AM = 3[1 + AW]E, 



|(p)AM, (37b) 
(37c) 
(37d) 



while the constraints ( 35b I and ( 35a I become respectively 



ff'-3E 



W'-3(p)AWr 
y 



y 



ff' - 3ET 



Y 

Y; 



0, 

= 0. 



(38a) 
(38b) 



The constraint (33b I leads to (38a) which is identical 
to ( 30 1 , so this constraint is automatically satisfied at all 
t (all 3 T) by a system of evolution equations based on 
the variables (p) and A' m -), 

Bearing in mind (17 1, (fl8|) and (26), the constraint 



( |33a[ ) in its form (|38bjcan be written as the integrability 
condition for I: 



Yi 



^ =r 



Further, since (39 1 holds and (36 1 implies 

-,1/2 



H =l = 



KC 

~3 



2 <P)-J( 3 ^) 



0. 

(39) 
(40) 



the constraint ( 38b ) implies (with the help of ( 23 1 , (24 1, 
Q, and pi): 



1 jKC 



; (p)AM - i( 3 ^)A( fc ) 



1/2 



(41) 



If the system (37) is self-consistent, then ( [41) must hold 
for every 3 T and so it must comply with ( 37b I . In- 



serting ( 41 1 into ( 37b I and eliminating time deriv ative s 



of (p), A"W, H, (*Wf and A< fc ) with the help of (37a 



(37c), (|37d| and the scaling laws ( 23 1— ( 24 1, it is straight 



forward to see that ( 37b I is identically satisfied at every 

3 r. 



DIMENSIONLESS VARIABLES: A 
DYNAMICAL SYSTEM 



In order to transform system (37a)-( 37d| into a dy- 
namical system, we need to recast H, E and (p) in terms 
of dimensionless "expansion normalized" variables. [3]. 
However, instead of using the expansion parameter in 
( 25 1 , it turns out to be easier to use H defined in d36| . 



This leads to: 



n 

(n) 
s 



H 
iV 

™ 2 (p) 

E 
H' 



(42) 

(43) 
(44) 



where Hq is a characteristic length scale. Notice that we 
can also define a dimensionless f2 parameter for the local 
density, which would be related to that of (43 1 by: 



KC 2 p 

3H 2 



(fi)[l 



v( m )l 



(45) 



We also need to make the "dot" derivative d/dct a di- 
mensionless operator, so we define: 



d_ 



J_d_ 



(46) 



Inserting (42|-(46l into ( 37a I— ( 37d I we obtain the fol- 



lowing dimensionless system: 

n, T = - [i + |(n>] n 2 , (47a) 

S, T = H[ S(3S - 1) + ± (A(™) + S) (ft)] , (47b) 
(0) >r = H{il) [<fl) -1], (47c) 
3HS [i + AW| . (47d) 



a(™) 



where jT = d/dr. 

The form of these equations suggests that further sim- 
plification follows if the three equations (47b I— (47dl can 



be decoupled from (47a 
nate transformation: 



Defining the following coordi- 
r = r(£,f), r = r (48) 



so that for any function A — A(t, r) we have A(t, r) — 
A(£(t, f),f) and all partial time drivatives in (47 1 can be 
expressed as 



'OA' 


OA 








9t_ 



= H 



dA 



where £ is selected so that: 
9i 1-1 



C = ln€, 



(49) 



(50) 
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The introduction of the new variable £ defined by ( 49 ) 



and ( 50 I removes the dependence on ft of the evolution 



equations <\47b} , fl47c[ ) and p7d l 

as; 

5£ 



5£ 



5 

3S 



(3S-l) + ^(A(" l ) + s) (0), (51a) 
1 + A (m) 



(n) [(O) - 1] 



while (47a I becomes 



c9£ 



(In ft) 



(51b) 
(51c) 

(52) 



The system (51) is formally a dynamical system con- 
structed with "expansion normalized" variables [3] . How- 
ever, it is still a system of PDE's, even if it only contains 
derivatives with respect to £ and r is basically a param- 
eter. As we show in Appendix B, this system of PDE's 
equivalent authonomous ODE system with restricted ini- 
tial conditions (restricted in order to fulfill the spacelike 
constraints) . 



VI. INITIAL CONDITIONS 



Bearing in mind Q, ([TTJ, (^-(pl) and (|40|-(glJ, 
intial conditions for the system ( |51[ ) can be constructed 
through the following steps: 

1. Choose the topological class of the initial hyper- 
surface 3 Ti in terms of the number of RSC (see 
Appendix A). A convenient choice is: 



(53) 



where f(r) is a (at least a C 2 ) function whose zeros 
correspond to RSC's and Hq is the same character- 
istic length scale as in ( 42 I and ( 46 1 . Strictly speak- 



ing, the choice of / is a choice of radial coordinate 
(see Appendix A). 

Construct dimcnsionless quantities out of pi and 

3 n, : . 



KC 2 Pi 

w 



3 n t 



(54) 



Obtain the orbit volume averages and contrast 
functions 



(m») 


~ 3# 2(Pl) ' 


(55a) 


(h) 


= 6^ <™> 


(55b) 


A ( m ) 

i 


mi 
(mi) 


(55c) 






(55d) 



4. The initial forms of the remaining variables are 
then: 



ft, = [(mi) - (h)] 1 ' 2 
(fli) = 



S, 



[mil 
(m^ - (ki) ' 

1 K) Ag - (h) Af 

2 ( mi ) - (h) 



(56a) 
(56b) 

(56c) 



Just as we proceeded with ( 37 1 , the fulfillment of con- 
straints ( 38a| and (38b I imply that the functional forms 
of ft, (ft) and S are the generalization of (56 1 for all times, 
that is: 



ft = 
(fi) = 

S = 



\(m)- (k)} 1/2 , 

( m ) 
(m) — (k) ' 

1 (to) A(™) - (k) AW 
~2 (m) 



where 



Hp) 



3H 2 



(k) 



(k) = 



(*n) (h) 



6H 2 



(57a) 
(57b) 

(57c) 



(58a) 
(58b) 



Inserting (57 1 and (58 1 into (51 1 and eliminating deriva- 
tives of (mYand (k) by means of p4|, dgg), ([37a]), fl37c 



and (37d), we can see that (57 1 and([58|) fully solve (47 1 



Thus (47) propagates the initial data given by ( 53 1— (56c 
with the constraints ( 38a I and ( |38b[ ) satisfied for all ^ T 



The standard approach to LTB dust solutions is based 
on solving the evolution equation pi). Therefore, it is il- 
lustrative to re-write such equation in terms of the initial 
value variables that we have defined here: 



e T = n 2 



(Cli 



(59) 



This form of ^ shows how the sign of (O^) — 1 determines 
the type of evolution of dust layers, just as the sign of 
K does it in d3|), leading to "parabolic" ((fij) — 1 = 0), 
hyperbolic" ((ft) - 1 < 0) and "elliptic" ((ft*) - 1 > 0) 
evolution like setting K = 0, K < 0, K > in See 
[51 |B] . Analytic solutions of ( 59 1 are given in Appendix 
C. 

In the numerical examination of dust configurations we 
will need to solve numerically system ( 37 ) , besides system 
( 51 1 . In terms of the variables and initial conditions ( 53 )- 



(58 1 the system (37 1 takes the dimcnsionless form: 



(ro), r 



-H 2 



(to) 
2 

3 s 2 - 2ft s - 



\m) 



-3<m)ft, 
3(1 + A( m ')s, 



(60a) 

(60b) 
(60c) 
(60d) 
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have been defined in ( 42 1 , (461, (58) Since (SI) > 0, this constraint (together with (63al) has 



where r, TL and 
and s = T,/Hq. This system will be specially useful for 
looking at the collapsing stage of re-collapsing configu- 
rations, a feature that cannot be studied by ( |5l"| because 
maximal expansion takes place as TL — * 0, and so (SI) di- 
verges and the solution curves cannot be extended further 
to study their collapsing stage. Also, (60 1 yields directly 
quantities like (to) and TL, while the metric functions I 



and r in (29 1 follow from the variables in (60 1 as: 



(to) 



1/3 



r = 



i 



i 

AM ' 



(61) 



and the so-called "curvature radius" Y = 
easily found as Y = Yi £. 



I gee can be 



VII. SOME QUALITATIVE AND ANALYTIC 
RESULTS. 

The solution curves [S(£, r), A< m )(£, r), (fi)(£, r)] of 
(51 1 evolve in a phase space which is a region of K 3 
parametrized by the coor dinates [S, A<- m \ (SI)}. The evo- 
lution variable is £, while each curve represents a funda- 
mental comoving observer labelled by a constant value 
r. It is highly desirable that the phase space coordinates 
remain bounded as these curves evolve along their maxi- 
mal range of £ for all r. We discuss in this section several 
useful analytic results and the relation between the phase 
space variables and their initial values. 

Considering that £ = In £ and 



(Si) - 1 = 



(h 



(m t ) - (ki) ' 
(k) 



(k) 



(62a) 
(62b) 



the functions TL, (SI) and S in (57l-(|58| can be given as 
(Sli) 



{n) ~ (Sli) - [(Sli) - l]eC' 
TL 2 = H?[(Sli)e- 3 t-((Sli)-l)e- 2 t] 

(ii,)AW-[(il,)-l]AWe^ 
2[(Sl i )-[(Sl i )-l]ei] ■ 



(63a) 
(63b) 
(63c) 



where we are emphasizing that these are now functions 
of (£, r). Notice that even if, in general, the surfaces of 
constant £ do not coincide with the 3 T (surfaces of con- 
stant t or t) , the initial hypersurface £ = does coincide 
with the initial Z T.;. 



A. Constraints on (SI) 



From(57l, we have the following important constraint: 

1 



/ON _ 1 <*> = _M- 

1 ; W 2 ^ 2 H 2 - ^ 2 7i 2 



(64) 



important qualitative consequences: it defines invariant 
subsets associated with (SI). If initial conditions for a 
range of r are selected so that (Sli) = 0, < (Sli) < 
1, (fii) = 1, or (Sli) > 1, then all solution curves for 
such range respectively comply with (SI) = 0, < (SI) < 
1, (SI) = 1, or (f2) > 1 for all their maximal extensibility 
range of of the evolution parameter £. 

Notice from (63 1 that, irrespectively of the sign of 
(Sli(r)) — 1, we have (SI) — > 1 as £ — * — oo (or, equivlently, 
as I — > 0). Thus, irrespectively of the type of dynamics 
i np0| given by the sign of (Sli) — 1, all solution curves of 
( |51 p will approach the plane (SI) = 1 near the initial big 
bang singularity (in the asymptotic range: £ — > — oo or 
*-»0). ' 

Since (fi) > 0, (64 1 implies that all solution curves 
with initial conditions < (Sli) < 1 will be constrained 
to the invariant subset given by the region < (SI) < 1 
for all values of £, while (63 1 implies that the such curves 
will evolve from the plane (SI) = 1 in £ — > — oo (I — > 0) 
towards the plane (SI) = in the asymptotic range £ — > 
oo (or £ — > oo). 



Equations ( 63 1 place no restriction in the range of £ 
of solution curves when (Sli) < 1, but curves with initial 
conditions given by (Sli) > 1 the evolution parameter £ 
will have a bounded range of maximal extendibility given 
by 



(SI) — > oo, as £ 



In 



(ft; 



(n 4 



i 



(65) 



Such solution curves start from the plane (SI) = 1 as 
£ — * — oo (i — * 0) and diverge for the limiting value of 
£ = ln^ given by (65 1. However, this is the same finite 



value of £ that makes H in (63b I vanish, so the blowing 
up (SI) — > oo of all solution curves in this case simply 
marks the values of £ associated with the "maximal" ex- 
pansion of dust layers (maximal value of £) in the "turn 
around" before the collapsing stage. Since (SI) blows up 
the collapsing stage cannot be described by ( |5l"j ) (though 
it can be described by ( 60 1 ) . 



B. Constraints on S and A (m) 

In order to appreciate the relation between the coordi- 
nates of the phase spa ce, \ S,A^ m \ (SI)], it is very useful 
to express 5* given by (57c) as: 

S-- 



i{{n>AW-[<fl)-i] aw} 



(66) 



where, from the scaling laws (23 1 and (24 1, we have 

1 -+- A (m ) 

A (m) = 1, (67a) 

= 2/3 + A > ] - I (67b) 
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Thus, for finite A^ and Af, it is evident that a suf- 
ficient condition for A^ m ' and A( fc ' to remain bounded 
for all £ is that initial conditions are selected so that a 
shell crossing singularity does not emerge (see Appendix 
A and references [2j El EJ IZ] ) ■ That is, we must have for 
all solution curves marked by (£, r) 

r(£,r) > 0, no shell crossing singularity, (68) 



(fc) 



where T has been defined in ( 18 1 . Notice, from ( 63c ) , ( 66 ) 



and (67), that even if (68 1 holds, S (like (f2)) diverges 
as TL — > (maximal expansion) for curves with (f2j) > 
1, though it remains bounded if (f2j) < 1. In order to 
overcome the blowing up of S and (f2) in the cases when 
(f2j) > 1 we have used the system (60 1 to examine the 



collapsing stage of the solution curves. 

It is important to remark that for initial conditions 
for which both (68 1 and < (fij) < 1 hold, the three 
coordinates of phase space, [S, A^ m \ (fl)], of all solution 
curves remain bounded and restricted to a finite region 
ofM 3 . 

Initial conditions that guarantee the fulfillment of ( |68[ ) 
were derived in terms of original variables in reference [7 
(see also |6]), in terms of initial value variables dis- 
cussed here in [2] (see Appendix A). 




FIG. 1: Phase space The figure shows all critical points and 
the invariant subspaces given by the surfaces (Q) = 1 and 
A (m) = -1 (Schwarzschild vacuum case). Particular cases, 
such as FLRW, parabolic evolution and the vacuum case in 
Lemaitre coordinates, are also invariant subspaces shown as 
lines going from a source to a sink. The FLRW sub-case 
contains the Einstein-de Sitter (EdS) universe (fi) = 1 (the 
critical point Cs) and the Minkowski-Milne space-time (M) 
(f2) = (the critical point Ca). 



VIII. CRITICAL POINTS AND PARTICULAR 
CASES 

In the general case in which none of the variables 
[S, A^ m \(fl}] is restricted to take any special constant 
value associated with particular cases, the coordinates 



and nature of critical points of (51 1 are: 



(SI) = 1: 

Ci S=l/2, A (m) =-1, source, (69a) 

C 2 5= -1/3, A(™) = -1, saddle, (69b) 

C 3 S = 0, A^ = o, saddle, (69c) 

(SI) = 0: 

C 4 -5=1/3, A (m) = -1, saddle, (70a) 

C 5 S = 0, A (m) arbitrary, sink, (70b) 

C 6 S = 0, A (rn) = 0, sink, (70c) 



The nature of some of the critical points of ( |51| changes 
when [S, A^ m \ (Si)} take certain specific particular values 
that correspond to various space-times that are particu- 
lar cases of dust LTB solutions. We examine these par- 
ticular cases and their critical points in the remaining of 
this section. The phase space and the location of critical 
points and location of all particular cases is depicted in 
figure [1} 



A. Parabolic (or "marginally bound") evolution 



For initial conditions (f2j) = 1 equations (63 1 imply 
(SI) = 1 for all £. This is an invariant set associated 



with the so-called "parabolic" solutions of ( 59 1 and cor- 
responding to (k) = ( 3 TZ) = K = (see 0161). A full 
closed analytic solution (compatible with (57l) is given 
by: 



A 



(m) 



A (ro) (£,r) = ^, (71a) 

[1 + A 4 (m) ] e 3 «/ 2 - A| m) 



Sfor) = -^AH Kir ). 



(71b) 
(71c) 



so that solution curves remain for all £ in the line S = 
_A(™)/2 that lies in the plane (fl) = 1 (see figure [l}. We 
have all variables fully determined. It is straightforward 
to verify that ( [71] ) are fully compatible and equivalent to 
analytic solutions given in terms of £(ct, r) (see [5] and 
Appendix C). 

Critical points on this case are: 



Ci 



source, 



sink, 



(72) 



The phase space evolution goes along the line [S, —2S, 1], 
from Ci to C3 (respectively, past and future global at- 
tractors). Since the critical point C3 is characterized by 
conditions of homogeneity A' m ) = S = 0, it corresponds 
to the Einstein-de Sitter universe (dust FLRW model 
with zero spacial curvature). 
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B. FLRW dust models 

If A^" 1 ) = 5 = with unrestricted (f2), we have an 
invariant set associated with the particular homogeneous 
and isotropic sub-case of dust FLRW models. The criti- 
cal points are in this case: 



sink, 



(73) 



Depending on the sign of (fij) — 1, the phase space evolu- 
tion of these models will take place in the line [0, 0, (fi)]. 
If (fli) — 1 < the evolution goes between the source 
C 3 and the sink C 6 (past and future attractors). If 
(O,-) — 1 > 0, it starts at C 3 and goes upwards to 
(O) — * oo towards maximal expansion. In the collaps- 
ing stage, the point C3 acts as a sink. If (f2) = 1, we 
have the Einstein-de Sitter universe and the evolution 
is constrained to the point C3 with coordinates [0,0,1]. 
Also, we can identify the sink C6, whose coordinates are 
[0,0,0], as Minkowski space-time in the Milne represen- 
tation. 



C. RSC's 

The solution curve corresponding to a RSC is also char- 
acterized by A(™) =5=0 with varying (f2), thus the 
phase space evolution of that particular curve will also 
be along the line [0,0, (Q)], just like the FLRW sub- 
case. This is not surprising because a RSC is a privileged 
isotropic observer. 

The fact that the RSC's have a separate evolution 
from the rest of the solution curves prevents the criti- 
cal point Ci (a past and future attractor associated with 
the initial and collapsing singularities) to be global in all 
cases. However, Ci is a global attractor for all "off cen- 
ter" curves in configurations with RSC and for all curves 
in configurations lacking RSC (wormhole topology and 
vacuum limit). 



D. Vacuum case 

It is known that the LTB metric ([lj can describe the 
Schwarzchild-Kruskal space-time in terms of the world 
lines of its radial geodesic test observers with 4-velocity 
u a = Sq. The equation of motion of these geodesies is 
identical to ^ with M = M , so that M' = = p 
and the Ricci tensor vanishes, while the constant Mo 
can be identified with the "Schwarzschild mass" (see 
13 El IE])- The function 3 TZ is the scalar 3-curvature of 
the hypersurfaces of simultaneity of these geodesic con- 
gruences and the average scalar curvature is related to 
their binding energy per unit rest mass of the test parti- 
cle Ek(r) = -K/2 = (-1/12) fT^ 2 = (-1/2)1^). 

This particular case can also be constructed by assum- 
ing a Dirac delta distribution for pi so that a nonzero (p^ 
follows from the integral definition (10 1. Also, without 



resorting to a Dirac delta, this vacuum case follows from 
the 3+1 formalism in equations ([34]) by setting p — 0, so 
a definition of a nonzero average density follows from the 
Weyl tensor component in |7]) and (27 1 as: 
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(P) 



(m)H* 



-M /Y* 



(74) 



Since m = p = but (m) > 0, as long as T is finite (see 
section X) we have the following sufficient condition to 
characterize the Schwarzchild-Kruskal vacuum: 



A (m) = -1 for all (£,r) 



(75) 



which defines an invariant set, since all solution curves 
with initial condition = —1 will necessarily be 

constrained to the plane A^" 1 ' = — 1, parametrized by 
[S, (fi)]. The phase space evolution and critical points 
depends on the sign of (Hi), which depends in turn on 
the binding energy of the radial geodesies. 

For zero binding energy (ki) = = (fij) — 1 (the 
so-called Lemaitre coordinates for Schwarzschild space- 
time [5]) the critical points are the global future and past 
attractors: 



Ci 



source, 



sink, 



(76) 



and phase space evolution is the line [S, — 1, 1], from the 
source Ci to the sink C2. 

For geodesies with negative binding energy: (ki) > 
or (fij) > 1 (the so-called Novikov coordinates [9]) the 
critical points are the same as those of the zero energy 
case, but C2 is now a saddle. The phase space evolution 
of the curves begin at the global past attractor (source) 
Ci, approach the saddle C2 and go upwards towards di- 
verging (f2) at the maximal expansion, terminating again 
at the global future attractor (sink) Ci. 

For geodesies with positive binding energy ({ki) < or 
< (Hi) < 1) the critical points are: 



Ci 

c 2 



source, 
saddle. 



C 5 

c 4 



sink, 
saddle. 



(77) 



In this case the curves begin at the global past attrac- 
tor (source) Ci, terminate at the global future attractor 
(sink) C5, approaching the saddles C2 and C4. This 
case is examined in subsection F of next section (see fig- 
ure [18 1. 

We examine the phase space evolution of solution 
curves of several general case dust LTB configurations 
in the following section. 



IX. LTB DUST CONFIGURATIONS 

A. Expanding low density universe 

A low density perpetually expanding dust configura- 
tion follows by selecting (m,)(0) < 1 and negative (ki) 
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for all r defined in an initial hypersurface 3 7^ with an 
open topology having a RSC. This leads to (f2j) < 1 or 
to "hyperbolic" dynamics. Such a configuration can be 
realized with the following initial value functions 



0.9 



nij 

kj 
Yi 



1 + tan 2 r ' 
tan 2 r 



H Q 1 tan r, 



(78) 



so that we have an initial density clump ( — 1 < A^ < 

0) but a curvature void (0 < < 2/3). It is easy 

to show that these initial conditions yield (fij) < 1 for 
all r, so that < (f2) < 1 for all £. These conditions 



also comply with the no-shell-crossing conditions ( 93 1 



so dust layers monotonously expand perpetually from an 
initial big bang. 

As shown by figure [2] the solution curves start at the 
past attractor or source Ci (big bang), evolve towards 
the saddle C 3 and terminate at the future attractor given 
by the line of sinks C5 in the plane (f2) = 0. The evo- 
lution of the RSC is constrained to C3 for all time. The 
function A( m ) remains negative for all £, thus we have 
density clumps at all 3 T. 




0.08 
0.12 S 



FIG. 2: Low density configuration. Dust layers expand 
from the past attractor (source) Ci, denoting the big bang 
singularity, approaches the saddle C3 and drops towards the 
future attractor given by the line of sinks C5. Notice that 
( 78 1 implies that the evolution of the RSC is constrained to 
C 3 for all time. 



( — 1 < < 0) evolves for all curves into density voids 

with A^" 1 ) > 0. In this subsection we provide a numeric 
example of dust configurations discussed in reference |10j 
(see also HUHa). 

An initial dust clump transforming into a void emerges 
by choosing the following initial value functions 



max 



m o 



1 + a 2 , tan 2 r ' 



Y = 



moo = 0.05, moi = 0.9, ao = 2.0, 

. , fcoo 

01 + l + /3 2 tan 2 r' 
koo = -5.0, feoi = -1.0, 00 = 0.5, 
Hq 1 tan r, 



(79) 



which comply with the no-shell-crossing conditions ( |93| . 
The evolution of this configuration is similar to the pre- 
vious one: dust layers also expand perpetually ("hy- 
perbolic" dynamics) from an initial big bang and the 
hypersurfaces 3 T have an open topology with a RSC. 
However, the contrast function, A^ m \ now passes from 
negative to positive as £ increases, indicating that ini- 
tial clumps transform into voids. This can be seen very 
clearly if we plot (see figure (J3j> ) the normalized den- 
sity profiles along the hypersurfaces 3 T from the function 
p/p c = m/m c = (m)(l+A ( - m ))/m c , where the subindex c 
denotes evaluation along the RSC marked by r = 0. The 
density as function of r and r can be readily found by 



solving numerically the system ( 60 ) for initial conditions 
( |79| . This transition from clumps to voids is depicted by 
figure [3j 



log 10 (p/p c ) 
0.55 




B. Low density void formation 

An interesting variation of a perpetually expanding 
configuration, similar to the one presented in subsection 
A and figure [2] is that where an initial density clump 



FIG. 3: Low density voids: density profiles The picture 
depicts the normalized radial density profile p/p c for different 
hypersurfaces 3 T marked by constant r. Notice how an initial 
clump at r = evolves into a voids profile for latter r. 



This transition from clumps to voids can also be seen 
from the form of the contrast density function A^" 1 ) given 
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in terms of £ and r, displayed in figure [4] which shows 
how A'™ 1 ' is initially negative but becomes positive for 
all solution curves with r constant (except the RSC at 
r = 0) as £ increases. 





FIG. 5: Low density voids The evolution is similar to that 
of Fig 1, except that A' m ' becomes positive as all solution 
curves drop to the line of sinks C5 . 



FIG. 4: Density contrast function The picture depicts the 
function A' m ' passing for all < r < ty/2 from negative to 
positive indicating that an initial clump evolves into a void 
profile. Notice how A (m) = at the RSC r = 0, and asymp- 
totically at r — *■ 7r/2 (or Y; — > 00). For solution curves with r 
close to the RSC, A*-" 1 ' remains close to zero but is neverthe- 
less positive. 



The evolution of this configuration in phase space is 
shown by figure [5j Solution curves start at the past at- 
tractor (source) Ci (big bang), approach the saddle C 3 
and evolve towards the future attractor (line of sinks) C5 
in the plane (A) = 0. The RSC evolves from the source 
C3 towards the sink C6. However, now the curves evolve 
in such a away that all curves (save the RSC) hit the 
line of sinks C5 at positive values of A( m ) (passage from 
clumps to voids). 



C. Re— collapsing high density universes 

Re-collapsing high density configurations follow by se- 
lecting (mj)(0) > 1 and (hi) > for all r, leading to 
(A,-) > 1 for all r, thus (A) > 1 for all £. Dust lay- 
ers expand from an initial big bang, reach a maximal 
expansion and then collapse (the so-called "elliptic" dy- 
namics [2J [5J [B]). These configurations are compatible 
with 3 % having either a closed topology with two RSC, 
an open topology with one RSC or a wormhole topology 
without RSC (see Appendix A). We examine the closed 
and open case separately below, and the wormhole case 
in section E. 

The case with spherical "closed" topology can be ob- 



tained with the following initial value functions: 



mi 
Y 



moi 



H>1 



mpo - toqi 
1 + sin 2 r 



1 + sin 2 r 



H Q 1 sin r, 



m 00 = 3.9, m i = 1.1 
fc 00 = 2.9, feci = 0.22 

(80) 



which comply with (91 1 at the turning point r = tt/2. 
The phase space evolution is illustrated in figure [6] The 
solution curves start at the past attractor (source) Ci 
(big bang) and evolve upwards as (A) — > 00, which cor- 
responds to maximal expansion given by TL — ► and 
£ -> ln[(Oi) /((ft,) - 1)]. The RSC also evolves upwards 
from the source C3. For off center solution curves this 
critical point is a saddle. 

Since the collapsing stage cannot be described with 
solution curves of system (51 1, we can examine 
this stage with system ( 60 1 by plotting the curves 
[S , (r,r),AM(r,r)^)(r,r)]7^where S and (SI) are de- 
fined by ( |44| and ( 43 1 . As shown by figure [7J the curves 
come downwards from infinite (A) and S to the critical 
point Ci, which is now a future attractor or sink, associ- 
ated with a second curvature singularity ( "big crunch" ) . 
As in the expanding state, curves near the RSC approach 
the saddle C3, while the RSC evolves along the line 
[0,0, (A)] from infinite (A) to (A) = 1 at the sink C 3 . 

The case with open topology with one RSC can be 
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FIG. 6: Phase space evolution of a high density re- 
collapsing universe with closed topology. Solution 
curves evolve from the past attractor or source Ci, approach 
the saddle C3 and diverge as £ — > ln[(fij)/((f2i) — 1)], mark- 
ing the maximal expansion of dust layers. Only the range 
< r < tt/2 is shown, with the RSC at r = evolving from 
the source C3 upwards. The off-center layers approaching 
the saddle C3 are those with r w while those further away 
from C3 are close to the "turning value" r — tt/2. Layers 
marked by tv/2 < r < n, including the second RSC at r = ir, 
would have identical evolution as those shown in the figure, 
with those marked by r ~ n near the saddle C3. 




FIG. 7: Phase space evolution of the collapsing stage 
in the re— collapsing universe with closed topology. 

Only the end stage of the collapse is described by plotting 
the phase space variables by means of system ( 60 1 . Solution 



curves evolve from infinite values of (fl) and S to the point 
Ci, which is now a future attractor or sink. The curves also 
approach the saddle C3. Only the range < r < 7r/2is 
shown, with the RSC at r = evolving from infinite values 
of (Q) to the C3, which is now a sink. Layers marked by 
tt/2 < r < n, including the second RSC at r = tt, would have 
identical evolution as those shown in the figure, with those 
marked by r w n near the saddle C3. 



obtained with the following initial value functions: 

5.0 

Y 



1 + tan 2 r ' 



2.0 



1 + tan 2 r ' 
H^ 1 tan r, 



(81) 

As in the case with closed topology, dust layers emerge 
from a big bang singularity reach a maximal expansion 
and then re-collapse in a big crunch, but now dust layers 
far from the RSC at r = re-collapse in very long time 
which becomes infinite in the limit r — *■ 7r/2. 

The evolution of the solution curves of system pi) in 
phase space is shown by figure [8] This evolution is similar 
to the closed spherical case depicted in figure [6] Since 
initial conditions yield (f2j) > 1 we have (f2) > 1 for all £. 
The curves (except the RSC) start at the past attractor 
or source Ci (big bang), approach the saddle C3 and 
diverge upwards as £ — * ln[(fij)/((f2,)— 1)], corresponding 
to maximal expansion (TL 0). The RSC evolves from 
the source C3 towards (f2) — * 00. 




FIG. 8: Phase space evolution of high density re- 
collapsing universe with open topology. Since (fi) — > 00 
as £ — > (flj)/((f2j) — 1), we plot arctan(f2), The curves start 
at the past attractor or source Ci (big bang) approach the 
saddle C3 and diverge upwards (maximal expansion). The 
RSC evolves from the source C3 to (f2) — > 00. The collapsing 
stage is not described. 



As in the case with spherical topology, the collapsing 
stage cannot described with pi) . We use then the solu- 
tion curves of ( 60 1 to plot the phase space variables in the 
collapsing stage in figure [9] which is very similar to figure 
[7] the critical point Ci is now a future attractor or sink 
associated with the collapsing singularity ( "big crunch" ) , 
while curves near the RSC approach the saddle C3 and 
the RSC evolves from infinite (f2) downwards to (f2) = 1 
along the line [0,0, (Q)]. 
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FIG. 9: Phase space evolution of the collapsing stage in 
the re— collapsing universe with open topology. Phase 
space variables are plotted by solving system ( 60 1 . Solution 



curves evolve from infinite values of (Q) and S to the point 
Ci, which is now a future attractor or sink. The RSC at 
r = evolves from infinite values of (Q) to the sink C 3 . 



with expanding curves filling the asymptotic range of r. 
For the collapsing curves the function TL obtained by 



solving the system ( 60 ) for initial conditions ( 82 I changes 



sign (passes from positive to negative). This function is 
plotted in figure [10] showing collapsing and perpetually 
expanding layers respectively in red and blue colors. 




0.0-5 



D. Structure formation scenario 

A very important type of configuration, associated 
with structure formation scenarios, follows by having a 
local region around a RSC evolving from a big bang sin- 
gularity and collapsing into a black hole (second singular- 
ity), while the region "outside" expands perpetually as a 
sort of "cosmic background" . This type of configurations 
have also been examined with the LTB original variables 
(see [TJ |5J [TSJ HH [HI US] and references quoted therein) . 

The structure formation scenario requires that (f2$) — 1 
be positive at the RSC and becomes negative for large r. 
It can be constructed by choosing > 1 near the RSC 
that becomes m, < 1 asymptotically, while ki must be 
positive at the RSC and become asymptotically negative. 
Such a choice will yield for re-collapsing curves (f2j) > 1 
near the RSC, so that (ft) < 1 and diverges for £ — > 
In[(f2j)/((f2j) — 1)]. For expanding curves we have < 
(ft,) < 1 and < (ft) < 1 for all £. The following initial 
value functions fulfill these conditions: 



ki 
Y, 



TO01 ~\ 

fcoi + 



m 00 - ropi 

1 + r 6 
fcpp - fcoi 
1 



r 6 



Ho 1 r, 



m oa = 1.5, m i = 0.2 
k a = 0.5, fern = -0.5 

(82) 



We have depicted in all the graphs the collapsing solu- 
tion curves in red and the perpetually expanding ones in 
blue. The collapsing curves evolve as the solution curves 
of figures [6] and [8] except that now these curves only 
cover a bounded range of r around the RSC at r = 0, 



FIG. 10: The function TL in a structure formation 
scenario. This function is obtained by solving the system 



(601 for initial conditions (82 1. For re-collapsing curves with 



(Qi) > 1 (red color) this function passes from positive (ex- 
panding) to negative (collapsing), while for perpetually ex- 
panding curves with (Qi) < 1 (blue) it remains always posi- 
tive. 



It is also interesting to visualize the the dimensionless 
density m (see figure 111. For collapsing curves (red) 
around the RSC at r 



this function decreases from in- 
finite values as dust layers expand from the initial singu- 
larity, reaches a minimal value (maximal expansion) and 
then increases again as curves go into a collapsing singu- 
larity. However, for perpetually expanding curves (blue) 
this function decreases monotonically from the initial big 
bang. 



The evolution of the solution curves of (51) in phase 



space is more complicated than for previously presented 
configurations. Collapsing curves (red) with (O^) > 1 
evolve as curves in re-collapsing high density configura- 
tions of figures [6] and [8j starting at the past attractor 
or source Ci (big bang), approaching the saddle C3 and 
diverging upwards ((Q) — > 00) as £ — > ln[(fij)/((fij) — 1)]. 
Perpetually expanding curves (blue) with < (Sli) < 1 
resemble the curves of the low density universe of fig- 
ure [2] starting at the same past attractor or source Ci 
(big bang) approaching the saddle C 3 and descending 
towards the future attractor given by line of sinks C5. 



This is shown in figure s [T2| and 13 

As shown by figure |12| all curves start at the past 
attractor Ci (big bang), those near the RSC (r = 0) ap- 
proach the saddle C3 (see figure 12 1 while those near the 



split from collapsing to expanding approach the saddle 
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FIG. 11: Density for a structure formation scenario. 

The initial big bang singularity is now associated to infinite 
density that decays as dust layers expand. Curves (red) with 
(Qi) > 1 reach a minimal density (at maximal expansion) 
and then density rebounds to infinity at the collapse. For the 
remaining curves with (Qi) < 1 (blue) density just dilutes 
monotonously as r grows. 




FIG. 12: Phase space picture of the structure forma- 
tion scenario. Red color corresponds to collapsing curves 
with (Qi) > 1, while perpetually expanding curves with 
< (Qi) < 1 are depicted in blue. All "off center" curves 
start at the past attractor Ci, perpetually expanding curves 
descend into the future attractor given by the line of sinks C3, 
collapsing curves approach the saddle C3 and go upwards to 
diverging (Q) (maximal expansion before collapse). The sad- 
dle C2 splits expanding from collapsing curves. The curve 
with (fij) = 1 (not shown) would go from Ci towards this 
saddle. 



C2 (see figure |l~3|). All collapsing curves end up diverging 
upwards (f2) — * 00 as £ — > ln[(f2j)/((fij) — 1)] (maximal 
expansion). Blue curves with with < (f2j) < 1 are per- 
petually expanding with < (Hi) < 1 and so < (Q) < 1 
for all £. Expanding curves also start at the source Ci 
(big bang), but do not approach the saddle C 3 , instead 
they approach the saddle C 2 and end up in the line of 
sinks C5. The RSC evolves from the source C3 upwards. 
Both the collapsing and expanding curves approach the 
saddle C2 , which splits the former curves going upwards 
to (Q) — > 00 and the latter downwards towards the line 
of sinks C5. The curve with (f2j) = 1 remains in the 
plane (f2) = 1, evolving from Ci towards C2- See figure 

on 

The collapsing curves behave near the collapse as the 
curves in figures [7] and [9] We show in figure [14] the col- 
lapsing stage of the curves plotted in figure [13] for a range 
of r around values where (f2j) — 1 passes from positive 
(collapsing curves) to negative (expanding curves). No- 
tice how collapsing curves end up falling into the critical 
point Ci, which is a future attractor or sink. As the 
curves approach the splitting saddle C2, the collapsing 
ones go to infinite (Q) and S, while expanding ones fall 
into the sinks C5. Because of the restricted domain of r 
that we used to illustrate the behavior near C2, the fig- 
ure doesn't show curves near the RSC approaching the 
saddle C3 and the RSC evolving from infinite (f2) into 
C3, which is a sink for this worldlinc. 




FIG. 13: Close up of the solution curves near the sad- 
dle C2.Red color corresponds to layers with (Qi) > 1 and blue 
ones to (Qi) < 1. The plot only depicts solution curves that 
correspond to r values near the change from re-collapsing to 
perpetually expanding behavior (where H becomes negative) . 
Notice how this saddle splits collapsing curves going upwards 
to (Q) — > 00 while expanding ones go downwards to line of 
sinks C 5 . Notice how the curve corresponding to (Qi) — 1 
remains in the plane (Q) — 1, starting from the source Ci and 
ending in the saddle C2. 
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FIG. 14: Collapsing stage for curves of figure |13| The 

curves come from solving ( 60 1 for initial conditions ( |82[ ). The 
end state of collapsing curves (red) is the point Ci, which 
is now a sink. Collapsing curves go upwards taking infinite 
values of (Q) and S. Notice how the saddle C2 splits col- 
lapsing curves going upwards to (fi) — > 00 while expanding 
ones (blue) go downwards to the future attractor given by the 
sinks Cs- 



E. Re— collapsing wormhole 



In dust configurations without a RSC ("wormhole" 
topology) , the 3 T can be homeomorphic to either § 2 xl 
or to § 2 x S . These configurations were examined in 
reference [T7] (see also [5J [5]) and can be constructed 
with Yi(r) having no zeroes in all its domain. Since any 
choice of such function will necessary fulfill Y/(r*) = 
for at least a "turning point" r = r* , these configurations 
can only be regular if they comply with the regularity 
condition ( |9"T| ), which implies a re-collapsing "elliptic" 
dynamics with (f2j) — 1 > (sec Appendix A). 



In the integral definitions of (p) and ( 3 1Z) in ( 10 1 and 
pT| we took the RSC as the lower bound of the integrals. 
This is a boundary condition on these definitions ensuring 
that regularity at that RSC is fulfilled [2EJE]. If there 
are no RSC then this boundary condition can be specified 
by demanding that rrij and fcj vanish at an asymptotic 
value of r along the 3 %. Such value can be then taken 
as the lower bound for the integrals. However, while 
this solves the mathematical problem, without a RSC 
the interpretation of (p) and ( 3 1Z) as average functions 
becomes less clear. 



A choice of initial value functions for the re-collapsing 



wormhole is given by: 

m 00 (l + a ) 



1 + ceo sec 2 r 
Ml+^p) 
l + 0o se ° 2 r ' 



moo = 1-1) ckq = 5.0 



too 



0.48, A) = 3.0 



Yi = m 1 



sec r, 



(83) 



Notice that the "turning point" (or "throat") in this ex- 
ample is located at r = r* = with (91 1 taking the form 
(ki)(0) — 1. Also, both to, — > and fc, — > as r — > ±7r/2. 
Therefore, we can define and (hi) by integrals like 



(10 1 and (111 by taking the lower integration limit at 



r = — 7r/2 and integrating forward in r up to n/2. The 
metric function Y = YJL obtained from ( 60 1 for this case 



is depicted in figure 15 




FIG. 15: Curvature radius of dust layers with worm- 
hole topology. Notice how for each hypersurface 3 T with r 
constant, Y — Yii has no zeros in its regular range and takes 
its minimal value at the "throat" at r = 0. Dust layers near 
the "throat" re-collapse in very a short interval of r, while 
layers further away re-collapse in much longer times. 



The evolution of the wormhole configuration in phase 
space is depicted in figure |16| Only curves in the range 
< r < 7r/2 are displayed, since (because of the sym- 
metric construction of the example) curves in the range 
— ir/2 < r < behave identically This evolution is sim- 
ilar to that of high density re collapsing universes, with 
all solution curves starting at the past attractor or source 
Ci associated with an initial big bang singularity and ris- 
ing towards maximal expansion (f2) — > 00. The lack of a 
RSC is evident and this means that Ci is a global past 
attractor. 

It is interesting to notice how the "throat" and curves 
near it rapidly leap into diverging (f2) because they reach 
maximal expansion and collapse very fast (see figure 15 1, 
but curves close to r = tt/2 spend a long time very near 
the plane (f2) = 1 and approach the saddle C3 where 
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FIG. 16: Dust wormhole There is no RSC. Notice how 
curves near the "throat" at r = rapidly reach maximal ex- 
pansion with ($7} diverging, while curves close to the asymp- 
totic value r — n/2 take a long time to reach this stage and 
approach the saddle point C3 with A^™' = 5 = 0. The lack 
of a RSC implies that Ci is now a global past attractor. 




FIG. 17: Dust wormhole: collapsing stage The figure 
describes the expanding and collapsing stages. The curves of 
figure [16] near the "throat" r — can be recognized and seen 
going upwards into infinite values of (fi) and 5, all in order 
to plunge downwards into Ci. Notice that the lack of a RSC 
makes this critical point a global past and future attractor. 



S = A 1 ™' = 0. We took as example a very simple sym- 
metric form of 3 % with only one turning point, but it 
is possible to choose any number of such points [T7]. It 
is also possible to construct a configuration whose hy- 
persurfaces 3 T are homeomorphic to torii § 2 x S 1 (see 

[51 E|). 

The collapsing stage for wormhole configurations is 
very similar to the stages depicted in figures [7] and [9] 
We show in figure [17] the expanding and collapsing stage 
for curves near the "throat" r = of this configuration. 
It is easy to recognize some of the curves plotted in figure 
16 which continue towards infinite values of (f2) and S, 
in order to plunge downwards to end in the critical point 
Ci, which is a future attractor or sink in this stage. 



F. Vacuum limit 

The evolution in phase space takes place in the plane 
A( m ) = — 1. As we commented in the previous section, 
for zero binding energy (f2j) = 1 (Lemaitre coordinates 
[5J) it goes from source Ci to sink C2, while for negative 
binding energy (f2j) > 1 (Novikov coordinates [§]) it goes 
from source Ci, approaches saddle C2 and goes upwards 
to (0) — > 00. The phase space evolution in this case is 
qualitatively analogous to that of a wormhole topology, 
which is not surprising since the Schwarzschild-Kruskal 
space-time has no RSC, hence Ci will be a global past 
and future attractor. 

We provide here an example with positive binding en- 
ergy (fli) < 1. As shown in figure [l8[ the solution curves 
emerge from the global past attractor source Ci, some 



curves approach the saddle C2 and some approach the 
saddle C4 and all terminate in the global future attractor 
or sink C5. 




FIG. 18: Vacuum case The solution curves correspond to 
the Schwarzschild-Kruskal space-time in coordinates given by 
radial geodesic congruence that is perpetually expanding with 
positive binding energy (thus (fii) < 1). Solution curves are 
confined to the plane A'*™' = — 1. The curves evolve from the 
global past attractor or source Ci towards the global future 
attractor or sink in C5 with 5 = 0, approaching the saddles 
C 2 and C4. 
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X. DISCUSSION AND CONCLUSION 



self similar sub-case of LTB metrics 



We have examined the class of dust LTB models in 
terms of the techniques gcncrically known as "dynamical 
systems" as they evolve in a 3-dimcnsional phase space 
parametrized by the variables [S, A^ m \ (fi)], whose phys- 
ical and geometric interpretation is intuitive and straight- 
forward. We can say in general that if initial condi- 
tions are selected so that "shell crossing" singularities are 
avoided (condition (68 1), then the phase space variables 
[S, A(™)] will be bounded. Also, the sign of (Qi)(r ) - 1 
for any given solution curve r = ro determines if it will 
evolve for all £ in the domain of J"o) — 1 that has 

the same sign. The various invariant subspaces, partic- 
ular cases and phase space evolution for representative 
dust configurations has been presented in detail in the 
previous sections. 

Except for the worldlines associated with RSC's, all 
"off center" solution curves begin at a past attractor 
given by a repelling source Ci, characterized by S = 1/2. 
(Q) = 1 and A^" 1 ) = —1 and associated with a big bang 
initial singularity. It is easy to show that Ci is indeed 
the only past attractor for all 'off center" solution curves: 
the asymptotic past regime for all curves is given by the 
limit £ — ► oo (or I — ► 0), so from (63a I we have (f2) — > 1 
in this limit. On the other hand, we can write A*™-* and 
AW as 



ds 2 



(m,-) r 



1, 



A (fe) = 



kj 



(ki)T 



- 1. 



(84) 



But, from the analytic forms given in Appendix C and 
for finite and regular m, and hi, we have in the limit 
I — > the following asymptotic behavior: T ~ l^l 2 for 
(fli) = 1 and T - l~ x for (pi) ^ 1. Thus, from @, 
we have A^" 1 ) — > —1 and A^ — > —1 in this limit, but 
then, from ( 66 ) and bearing in mind that (f2) — ► 1, we 
have S -►1/37 

As shown in figures [7] [9] [14] and [T7] all re-collapsing 
curves also terminate in Ci, which is in this case a fu- 
ture attractor for these curves and is associated with the 
collapsing ("big crunch") singularity. It is only the fact 
that RSC's do not evolve nor terminate (collapse) in Ci 
what prevents this attractor to be a global one, though 
in configurations lacking a RSC (wormhole topology and 
vacuum limits) this attractor is indeed global. 

In order to prove that near Ci we have a self-similar 
regime we use the definitions Q6J) , (17 1, (36 1 and (44 1. 
The condition S = 1/2 can be reduced to 



Y Cx = r(l + tf) 2 / 3 , 
r 



Cl (l + z?) 1 ^' 



with 



(85) 



where t = H$ct and we have eliminated two functions of 
the form a(ct) and b(r) by trivially re-labeling the time 
and radial coordinates. Since (Q,) = 1 implies (fcj) = 
= ( 3 7Zi), the metric (JT|) near Ci takes the form of the 



with 



-dr 2 + A 2 (d) dr 2 + r 2 B 2 ($) (dO 2 + sin 2 Odcf), 
A(0)=Y^, rB(0)=y Cl . (86) 



The fact that A( m ) = — 1 doesn't imply that Ci is a 
Schwarzschild vacuum, but that T — ► oo as t — ► 0. Eval- 



uating r for the self-similar metric forms (85) and (86 1 
with Yi = Y Cl (t i: r), we get 



r 



3 + ?9 l + 



3 + tf, l + d 



(87) 



with <&i — Ti/r, which shows how T — > oo at the locus of 
the big bang singularity $ = — 1 where Yci = rB(d) — > 0. 
Notice that the self-similar case has no shell crossing 
singularity since T vanishes at 1? = — 3 but the evolution 
range is d > — 1 (for collapsing stages we would have Y 
and Y' given by the same expressions as (851, but with a 



minus sign, so that the "big crunch" happens before the 
locus of T = 0). 

Therefore, we can consider the source Ci as marking a 
self similar limit and since all "off center" solution curves 
start at this critical point (associated with a big bang 
singularity), we can say that dust LTB models have a 
self-similar behavior near this singularity. Since collaps- 
ing curves in re-collapsing configurations also end in Ci 
(a sink) , we have self similar behavior near the collapsing 
singularity. 

Another interesting feature worth remarking is the role 
of the saddle C2 in the "structure formation" configura- 
tions discussed in section IX-D. As shown by figures [L3| 
and 14 this saddle splits collapsing curves (Oj) > 1 from 
the expanding ones < (f2j) < 1. The former evolve 
back into Ci (now a sink) while the latter fall into the 
line of sinks C5. This behavior denotes a basic instabil- 
ity of the invariant subspace of the "parabolic" evolution 
(ft) = 1, since for any solution curve in this space any 
arbitrarily small perturbation i5 < 1 on initial condi- 
tions (fij) = 1 + 6 will trigger a radically different "repul- 
sive" evolution away from (f2) = 1, either to a collapsing 
regime (6 > 0) or to perpetual expansion (5 < 0). 

The splitting saddle C 2 is marked by = -1, S = 

— 1/3 and (f2) = 1, it can be shown from the ana- 
lytic solutions in Appendix C and from numerical tests 
that r is finite as the curves approach these values (r 
diverges for later times as the collapsing curves actu- 
ally collapse). Therefore, the fact that A^ 1 ™) — * —1 
for curves approaching this saddle means that these 
curves experience near C2 unstable conditions similar to 
those of a Schwarzschild-Kruskal vacuum parametrized 
by Lemaitre coordinates made with radial geodesies with 
zero binding energy ((fc) = or (f2) = 1). In fact, the 
point C2 represents the global future attractor (sink) for 
this marginally bound vacuum configuration. 

The association of the saddle C2 with a Schwarzschild- 
Kruskal vacuum roughly conforms with an intuitive pic- 
ture. A very rough, hand waving, description of condi- 
tions near C2 could be: collapsing dust layers "retreat" 
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inwards towards a smaller collapsing region, while ex- 
panding layers "retreat" outwards in a fast expansion, 
thus we have a sort of unstable "evacuation" effect that is 
roughly similar to a Schwarzschild-Kruskal vacuum with 
the Schwarzschild mass being the accumulated effective 
mass Mq = (\/2)(rrii)Y 3 of the collapsing layers. It is 
known [5J [5] that in this type of configuration the more 
"external" layers in the collapsing region ((Hi) > 1) take 
infinite time to collapse, so that the resulting black hole 
perpetually takes accretion from these "border" dust lay- 
ers. However, we can assume that the mass contribution 
from this accretion is small and lengthy to come, so that 
once most "inner" layers with (f2j) > 1 have collapsed 
the expanding observers close to (f2,) = 1 do perceive a 
sort of approximate Schwarzschild-Kruskal vacuum with 
zero binding energy. 

We feel that the qualitative numerical treatment pre- 
sented in this article can serve for studying applications 
of astrophysical interest of dust LTB solutions. Perhaps 



for this purpose the system ( 60 1 can be more useful than 
the dynamical system (51 1. It is also interesting to see 



this methodology gcnerically applied to generalizations of 
these solution, such as LTB-de Sitter (dust with cosmo- 
logical constant) or the most general source compatible 
with the LTB metric: a mixture of dust and an inho- 
mogeneous and anisotropic fluid, which can respectively 
model cold dark matter and dark energy. Extensions of 
this work along these lines are presently under consider- 
ation. 
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APPENDIX A: Regularity and singularities in dust 
LTB solutions 

An extensive discussion of regularity conditions and 
important geometric features of dust LTB solutions can 
be found in [3 |6]. These features were also discussed 
with the variables of this article in [2] . We present here 
a brief review. 



Regular Symmetry Centers (RSC) and topology of 

the 3 % 



LTB dust solutions admit up to two RSC's, which are 
the regular wordlines made up by the space-time evolu- 
tion of the fixed point of the rotation group SO(3) where 
the orbits have zero area. In terms of the initial value 
functions, the RSC is the comoving worldline marked by 
a zero of the function Yi(r). Since the freedom to choose 
a radial coordinate in the LTB metric means that this 
function can be arbitrarily chosen, it is useful to select it 



as Yi — Hq 1 f(r), where f(r) is a (at least) C 2 real func- 
tion that describes in a simple manner the topological 
configuration or homeomorphic class for an initial 3 7l. 
We have the following choices: 

• "Open topology" : 3 % homeomorphic to M 3 , 



/(r)=tanr < r < tt/2 One RSC at r = 0. 



• "Closed topology" : 3 % homeomorphic to S 3 



/(r)=sinr < r < tt Two RSCs at r = 0, tt. (89) 



'Wormhole" : 3 % homeomorphic to S 2 x 



/(r) = secr - tt/2 < r < tt/2 Zero RSCs. (90) 



For the open topology f(r) is a monotonously increasing 
function. In fact, we could have just selected / = r 
for simplicity, but the choice / = tanr allows one to 
examine the asymptotic behavior associated with Y t — > 



oo in solution curves of (51 1 



In the wormhole topology the function Yi must be pos- 
itive and without zeroes. Numerical tests show that 
monotonous functions without (at least) one "turning 
value" (zero of Y() yield finite proper radial distances 
along the 3 % as r — > ±oo (depending on the choice of 
Yi). This indicates that such functions are inadequate 
choices of the radial coordinate, as they do not provide 
a full analytic extension of the manifold. Therefore, we 
select Yi > with at least a turning value r — r* within 
the open range of r. In the choice (90 1 this value is r* = 
(the "throat" of the wormhole) . 

In the "closed" topology we must have a turning value, 
which with the choice p9| ) is r* = tt/2. Since in both 
cases (closed and wormhole topologies) we have Y-(r*) = 
0, regularity of the metric function g rr 



Y' 



IVY' 



(l-JiQi/a [l-Hl(h)Y^ 

requires (see [H [51 [5] and [T!5]) that 3 lZi be selected so 
that K(r*) — 1, or equivalently: 



( 3 Ki)(r*)=H$(ki)(r*) 



(91) 



Given a choice of homeomorphic class for the 3 % all reg- 
ular 3 T have the same class (see [2]). All configurations 
examined numerically corresponding to closed or worm- 
hole topologies comply with (|9Xl) . 
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Singularities 



APPENDIX B: An equivalent dynamical system 



LTB dust solutions present two types of curvature sin- 
gularities marked by (ct,r) values such that: 

£(ct,r)=0, "bang" or "crush", (92a) 
T(ct,r)=0, "shell crossing" (92b) 

The "bang" or "crush" singularities are an inherent fea- 
ture and cannot be avoided, though it is always possible 
to select initial conditions so that "shell crossing" do not 
arise (either for all the evolution time or for a given range 

t > tO- 
Iiiitial conditions to avoid unphysical shell crossing sin- 
gularities were given in terms of the original LTB vari- 
ables by Hellaby and Lake [7] (see [5j [3]), and in terms 
of the initial value functions defined in this article by 
0. For broadly generic configurations that exclude very 
special (and unstable) features, we have: 

for < (Qi) < 1 : 

-1 < A| m) < 0, -2/3 < A| fc) < 2/3, (93a) 

for {Qi) > 1 : 

-1 < A| m) < 0, -2/3 < A| fe) < 0, 



We have examined the evolution of solution curves of 



2TTP, 

"3 



A ( m ) _ 

i 

P t Q t - 1 



3 A (fc) 
2 1 



> [PiQi - 1] A 



(m) 



Af ) > 0, . 



(93b) 



where 



Pi = 



x 



3-3/2 ' 

arccos(l — x) — \J x(2 — x), 
( 3 TZi) _ 2jfei) _ 2[(Qj) - 1] 

KC 2 p l 



{III, 



(Sit 



Equations ( |93| clearly show how initial conditions 
characterized by density and 3-curvature clumps (see 
( 15 1) lead generically to avoidance of a shell crossing sin- 
gularity if < (f2j) < 1, and as a corollary, to finite 
and bounded values of A'" 1 ' and S along solution curves. 
This regularity is more stringent for the case (f2j) > 1. 



However, even if ( 93 1 hold at the 3 7i it is still possible to 



have density or 3-curvature voids forming along the evo- 
lution of the curves, all without violating (68). This was 



discussed by [TU], but their methodology is too compli- 
cated. It is more clearly seen here from the fact that the 
scaling laws (23 1 and (24 1 do not, necessarily, forbid A^ 1 ™) 
and/or A( fe ) from reversing the sign of the initial A^ 

(k) 

and A^ . The configuration discussed in section IX-B is 
a numerical example of this type of "clump turning into 
void" configuration. 



the autonomous system (51 1 in the 3-dimcnsional phase 
space parametrized by [S~~A( m \ (fi)]. However, even if 
(51 1 has only derivatives with respect to the evolution 
parameter £, the differential equations are still PDE's. 
We provide below formal proof that such a system of 
autonomous evolution PDE's is equivalent to a system of 
ODE's with initial conditions restricted by the fulfillment 
of space-like constraints. 

Let JFjV) C M and J^te-r) C R be, respectively, the reg- 
ularity range of r and the maximal range of extendibility 
of £ for a given r, the system (51 1 can be associated with 
the flow <i>£ : R 3 — ► R 3 of the vector field d/d^ given by 



(95) 



where Xf : J 7 /^ — > R 3 is the initial state represented as a 
curve in R 3 given parametrically as 



X i (r) = [S i (r),A ( r\r),{Q i (r))}, 



(96) 



is a solution surface as- 



and X : !F(^ r ) X Fu.) — > 
sociated with X; by means of (|95|). This surface can be 



represented parametrically as the following surface in R : 

X(£,r) = [S{i,r), AW((,r), (fi(£,r)>], (97) 

so that Xi{r) = X(0, r). Given a point Xj(r*o) € R 3 there 
is a unique orbit or solution curve C ro : Ft&ro) ~ > R 3 of 
(51 1 which can be represented parametrically as: 



C ro (0 = X((;,r Q ), 



(98) 



so that C ro (0) = Xi(ro). The set of solution curves C T (£) 
for all r € generates the solution surface X(£,r). 
In order to relate (51 1 with a system of ordinary dif- 



ferential equations, we notice that every solution curve 
( 98 1 associated with a fixed r = r is a solution of the 



following system: 



dSo 



S (3S -l) + i(A(" l) + S ) (Q) , 



dA 



(m) 



— 3 So 



1 + A > 



(m) 



(Q) [(Q) - 1] , 



where the functions (So, Aq™\ (O,)o) '■ F(t,ro) 
given by 



So — 
A (m) - 



A (m) (C,ro) 



(99a) 
(99b) 

(99c) 
> R arc 

(100a) 

(100b) 
(100c) 
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and comply with the initial conditions: 

[5 (0), A { o m) (0), (fi)o(O)] = Xi(ro). (101) 
The fact that every solution curve of d5ll will be a solu- 



tion of a system like (99) for initial conditions given as 



in (101 1, suggests considering a system of ordinary dif- 



ferential equations associated with the same flow as ( 95 ) 



ds 
djj 



s(3s-l) + ~(6 + s) w, 



= 3 s [ 1 + 6] , 

for an initial state: 

f =[*(0),*(0),w(0)] = [«o,$o,wo] 



(102a) 
(102b) 
(102c) 

(103) 



and the solution surface given by a parametric form sim- 
ilar to (|97|: 



x= [s(C,s ), 5(£,6 ), u}(£,u )], 



so that x = x(0). Since the solution surface of this 
system can depend smoothly on initial conditions x*o, 



each one of the functions (s, 5, to) solving (1021 is a 



one parameter family of functions. In principle, there 



are many ways in which the initial state (103) can be 
parametrized. In particular, this state can be given 
as any curve [so(a), 6o(a), u>o(a)] in M 3 that intersects 
the solution curves x(£,x*o) only in one point, hence 
we can always find suitable one-parameter functions 
So (a) , <Jo (a) , u>o (a) such that x is expressible as: 

x(t a) = [«(£, ao(a)), 6 (a)), wo(a))]. (105) 



If, in particular, the initial state (103) is parametrized 



as: 



x (a) = [s (a),6 (a),wo(a)} = Xi(a), 



(106) 



with a £ J'(r), then for every a there is a C°° inclusion 
map ^ : K 3 — > K 3 mapping every solution curve of (51 ) 



given by ( 98 ) for a fixed r G ^"( r ) to a solution curve of 
the system 
is 



02 1 given by ( 105 ) with a fixed a = r. That 



V>(X(£,r))=x(£,a), 



and so, the solution surfaces of (51 1 are subsets of the so- 



lution surfaces of ( 102 ) obtained by restricting the initial 



states by the constraint (106 1 



Notice that ( 102 ) admits many solution curves that 



do not comply with (106 1, and so are incompatible with 
the constraints ( 53 1— ( 56c I and so they are incompatible 
with the evolution of an inhomogeneous dust source (LTB 
solution) complying with Einstein's field equations ex- 
pressed in terms of the 3+1 decomposition [3] under the 



regularity assumptions that we have made (see Appendix 
A). 



However, the solution curves of (1021 satisfying (106 



are equivalent (under ( 107 1 ) to solution curves of (51 
associated with the flow (951. Since (102) is a system 
of autonomous ODE's, geometric features such as criti- 
cal points and invariant spaces are well defined, so the 



qualitative analysis of (51 1 can be performed on (102) as 
long as we only consider its solution surfaces satisfying 
( 106 1. In other words, we can say that the dynamical sys- 



tem associated to inhomogeneous dust LTB solutions is 
a sort of "reduced" version of ( 102 1 with its initial states 
restricted by (106 1. Since, as we have shown in sections 



IV and V, the space-like constraints are satisfied at all 
space-like slices once they hold at the initial z %, the in- 
corporation of the radial dependence of the solutions only 
requires this restriction on the initial states. 

In order to simplify the notation we found it useful to 
keep the same names for the variables S,A^ m \ (fl) and 
to keep referring to ( |51[ ), with the understanding that 
the system we are really considering is ( |102[), t hus every 
mention of solution curves and surfaces of (15 ID refers to 



solution curves and surfaces of (102) mapped by (1071 



(104) and complying with (106) 



Appendix C: Analytic solutions 

Analytic solutions for LTB dust solutions follow from 
solving the evolution equation Q: 

Y 

Using the variables defined in sections III and IV, this 



equation becomes equation ( 59 1 

'<fii> 





2 r 




9t. 



m) 1) 



where: 



I 

2M 



Y 



T = Hn Ct 



(mi) Hq Y?, K^(k t )H^Y t 
(mi) (rrii) 



2 V 2 

i 1 



n 2 

(h) 



(rrii) - (hi 
(h) 



(107) 

We summarize the solutions of this equation below 



Parabolic solutions 



We have: (fti) = 1, so that (h) = and H 2 — (rrii). 
These are the only ones that can be given explicitly as 
I = £(t, r) from: 

2(^ 3 / 2 -l) 



T = Tj ± 



3Ti.i 



(108) 



21 



where the ± sign describes expanding (+) and collaps- 
ing (— ) dust layers, respectively characterized by r > 
T bb, t,T > an d t < r bb , £ T < 0. Notice that the "big 
bang time" r bb = Hoct hb is no longer an independent 
function but must be found by setting £ = in ( 108 1: 



(109) 



where now expanding and collapsing layers respectively 
correspond to "— " and "+" . 



Hyperbolic solutions 

We have: < (fi,-) < 1, so that (fe,) < and 
H\ = (rrii) + \(ki)\. In this case we cannot obtain 
I = £(r, r) in closed explicit form, but in parametric form: 
[£{t], r), t(j], r)] or implicitly as r — t(£, r). The solutions 
are given in parametric form as: 



£(r),r) 
r(T],r) 



(fli) [cosh 7] — 1] 
2(1 -(fli)) ' 
(tti) [(sinhr? - rj) - (sinhr?; - %)] 



(110a) 



2W 4 (1 - <^)) 3 / 2 



(110b) 



where 



T)i — arccosh 



((li 



2^1-^ ^ 
1 ) , sinh = v (.111) 



and ± denote expanding (+) and collapsing (-) layers. 
The coordinate locus of the central singularity (the "big 
bang time" function r bb ) can be obtained by setting rj = 
in Qll0b||, ( so that £ = 0) leading to 



r bb = n T 



(Oj) (sinh?7, - ?7;) 

2W, (1 - (fii))3/2 ^ 



(112) 



where now "-" and "+" respectively correspond to ex- 
panding and collapsing layers. 



Elliptic solutions 

We have: (fij) > 1, so that (k,) > and = 
(hi). The parametric solutions are 

(fli) [1 - cos rj\ 



£(r?,r) 
r(?7,r) 



2((n 4 )-i) ' 

(fii) [(?7 - sin 77) -(?/,- sin 774)] 



(113a) 



2H i ((n i )-l) 3 / 2 



(113b) 



where 



2y/Wi 



m = arccos _ !j , sin^= R) .(114) 

As with the hyperbolic case, the locus of the central sin- 
gularity (the " big bang time", r bb ), follows by setting 
77 = in (113b I, we get 



(fli) (jjj - sinrji) 
2H 4 «a)-l) 3/2 ' 



(115) 



Dust layers in elliptic solutions reach a maximal expan- 
sion as £. jT = in ( 59 1 , corresponding to 



(fi, 



(n< 



i' 

7T (f2j 



7/ = 7T, 



2w i «n i )- 1)3/2 



(116a) 
(116b) 



For £ :T < dust layers collapse a nd so a collapsing 
phase is given by parametric solution ( 1 13b I for ir < rj < 



2tt). The coordinate locus of the collapsing singularity 
("big crunch"), t 6c , follows from setting 77 = 2ir in (113b I: 



7T (f2i 



Wi^)- 1)3/2 



(117) 



Since the expanding and collapsing singularities are 
respectively marked by 77 = and 77 = 2tt, then the layers' 
evolution is symmetric with respect to 77 = n marking the 
maximal expansion (116 1. However, as shown by ( 115[ ), 
(116 1 and (1171, the maximal expansion time depends 
on the IVF's (f2j) and Hi, so it does not coincide (in 
general) with any 3 T (hypersurface of constant r) and 
the evolution for general IVF's (f2j), Hi is not "time- 
symmetric" with respect to the maximal expansion. 

While for the hyperbolic and elliptic cases we do not 
have an explicit closed form solution £ = £(t, r) . it is 
possible to use the parametric solutions presented so far 
to visualize graphically the function £(t, r) . This can be 
done by plotting the 3-dimcnsional parametric surfaces: 



[r(v,r),r,£(7],r)], 



(118) 



where the radial dependence is completely specified by 
the IVF's that enter in (TTlOh— dTTXl and p!3l-([TT4l. 
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